MOG RTX meta analysis
Bo Sun
2022-03-30
library(meta)
library(metafor)
library(ggrepel)
library(dmetar)
library(smd)
library(tidyverse)
library(PerformanceAnalytics)
library(ggplot2); theme_set(theme_bw(base_size = 15)+
theme(axis.text.x = element_text(angle = 45, hjust=1),
strip.text = element_text(colour = 'black', face="bold",size=12),
panel.grid.major = element_blank(),
panel.grid.minor = element_blank(),
panel.border = element_rect(size = 0.7),
axis.ticks.length=unit(.10, "cm"),
axis.ticks = element_line(size=0.7),
strip.background = element_blank()))##Loading data and prepping tables It is worth noting that using MD results in higher heterogeneity values than using SMD. Here I use MD as we have discussed. The significance does not change, but it affects the multi-model inference later. Using SMD, there is no significant influence… using MD there is.
Assessing bias
Here we see that MOG papers are biased. We will have to state this as a limitation
# Generate funnel plot (we do not include study labels here)
funnel.meta(mogmetagen)
title("MOG ARR")
# Plot over points
eggers.test(mogmetagen)## Eggers' test of the intercept
## =============================
##
## intercept 95% CI t p
## -2.654 -3.62 - -1.69 -5.404 5.847195e-05
##
## Eggers' test indicates the presence of funnel plot asymmetry.NMOSD papers are much more even.
# Generate funnel plot (we do not include study labels here)
funnel.meta(nmosdmetagen)
title("NMOSD ARR")
# Plot over points
eggers.test(nmosdmetagen)## Warning in metabias.meta(x, k.min = 3, method = "linreg"): 1 observation(s)
## dropped due to missing values## Eggers' test of the intercept
## =============================
##
## intercept 95% CI t p
## -0.591 -2.2 - 1.02 -0.72 0.4797325
##
## Eggers' test does not indicate the presence of funnel plot asymmetry.MOG EDSS scores are not biased across studies.
# Generate funnel plot (we do not include study labels here)
funnel.meta(mogedssmetagen)
title("MOG EDSS")
# Plot over points
eggers.test(mogedssmetagen)## Warning in metabias.meta(x, k.min = 3, method = "linreg"): 1 observation(s)
## dropped due to missing values## Eggers' test of the intercept
## =============================
##
## intercept 95% CI t p
## -0.265 -2 - 1.47 -0.298 0.7730529
##
## Eggers' test does not indicate the presence of funnel plot asymmetry.NMOSD EDSS scores are also not biased.
# Generate funnel plot (we do not include study labels here)
funnel.meta(nmosdedssmetagen)
title("NMOSD EDSS")
# Plot over points
eggers.test(nmosdedssmetagen)## Warning in metabias.meta(x, k.min = 3, method = "linreg"): 1 observation(s)
## dropped due to missing values## Eggers' test of the intercept
## =============================
##
## intercept 95% CI t p
## -0.257 -1.48 - 0.97 -0.411 0.6878221
##
## Eggers' test does not indicate the presence of funnel plot asymmetry.Forest plot MOG ARR
Here a random effects model reveals significant difference in mean ARR with lower ARRs post RTX.
##plot mog ARR forest plot
forest.meta(mogmetagen,
sortvar = yi,
predict = TRUE,
print.tau2 = FALSE,
leftlabs = c("Study", "g", "SE"))
print(mogmetagen)## Review: Mogmeta
##
## Number of studies combined: k = 18
##
## MD 95%-CI z p-value
## Random effects model -0.9213 [-1.2424; -0.6002] -5.62 < 0.0001
##
## Quantifying heterogeneity:
## tau^2 = 0.2650 [0.0813; 0.9016]; tau = 0.5148 [0.2851; 0.9495]
## I^2 = 87.2% [81.3%; 91.3%]; H = 2.80 [2.31; 3.39]
##
## Test of heterogeneity:
## Q d.f. p-value
## 133.01 17 < 0.0001
##
## Details on meta-analytical method:
## - Inverse variance method
## - Restricted maximum-likelihood estimator for tau^2
## - Q-profile method for confidence interval of tau^2 and tauForest plot NMOSD ARR
As shown before in Damata et al., NMOSD post RTX ARRs are significantly lower.
##plot mog ARR forest plot
forest.meta(nmosdmetagen,
sortvar = yi,
predict = TRUE,
print.tau2 = FALSE,
leftlabs = c("Study", "g", "SE"))
print(nmosdmetagen)## Review: NMOSDmeta
##
## Number of studies combined: k = 23
##
## MD 95%-CI z p-value
## Random effects model -1.7324 [-2.1507; -1.3141] -8.12 < 0.0001
##
## Quantifying heterogeneity:
## tau^2 = 0.7079 [0.3618; 2.6610]; tau = 0.8414 [0.6015; 1.6313]
## I^2 = 79.0% [69.1%; 85.7%]; H = 2.18 [1.80; 2.65]
##
## Test of heterogeneity:
## Q d.f. p-value
## 104.64 22 < 0.0001
##
## Details on meta-analytical method:
## - Inverse variance method
## - Restricted maximum-likelihood estimator for tau^2
## - Q-profile method for confidence interval of tau^2 and tauMerged ARR analysis
The merge analysis reveals that NMOSD patients benefit more than MOG patients. Between groups Q=9.09 p=0.003.
##plot mog and nmosd together
mergemetagen <- update.meta(mergemetagen,
subgroup = Group,
tau.common = FALSE)
forest.meta(mergemetagen,
sortvar = yi,
predict = TRUE,
print.tau2 = FALSE,
subgroup = TRUE,
leftlabs = c("Study", "MD", "SE"))
print(mergemetagen)## Review: mergemeta
##
## Number of studies combined: k = 41
##
## MD 95%-CI z p-value
## Random effects model -1.3866 [-1.6791; -1.0942] -9.29 < 0.0001
##
## Quantifying heterogeneity:
## tau^2 = 0.6122 [0.3644; 1.4990]; tau = 0.7824 [0.6037; 1.2243]
## I^2 = 93.2% [91.7%; 94.5%]; H = 3.85 [3.47; 4.27]
##
## Test of heterogeneity:
## Q d.f. p-value
## 592.56 40 < 0.0001
##
## Results for subgroups (random effects model):
## k MD 95%-CI tau^2 tau Q I^2
## Group = MOG 18 -0.9213 [-1.2424; -0.6002] 0.2650 0.5148 133.01 87.2%
## Group = NMOSD 23 -1.7324 [-2.1507; -1.3141] 0.7079 0.8414 104.64 79.0%
##
## Test for subgroup differences (random effects model):
## Q d.f. p-value
## Between groups 9.09 1 0.0026
##
## Details on meta-analytical method:
## - Inverse variance method
## - Restricted maximum-likelihood estimator for tau^2
## - Q-profile method for confidence interval of tau^2 and tauForest plot MOG EDSS
EDSS appeared to be significantly better in those post RTX.
##plot MOG EDSS
forest.meta(mogedssmetagen,
sortvar = yi,
predict = TRUE,
print.tau2 = FALSE,
leftlabs = c("Study", "g", "SE"))
print(mogedssmetagen)## Review: MogEDSSmeta
##
## Number of studies combined: k = 10
##
## MD 95%-CI z p-value
## Random effects model -0.8384 [-1.4146; -0.2621] -2.85 0.0044
##
## Quantifying heterogeneity:
## tau^2 = 0.4549 [0.0231; 3.3570]; tau = 0.6744 [0.1520; 1.8322]
## I^2 = 55.4% [9.4%; 78.1%]; H = 1.50 [1.05; 2.14]
##
## Test of heterogeneity:
## Q d.f. p-value
## 20.20 9 0.0167
##
## Details on meta-analytical method:
## - Inverse variance method
## - Restricted maximum-likelihood estimator for tau^2
## - Q-profile method for confidence interval of tau^2 and tauForest plot NMOSD EDSS
EDSS is significantly improved post RTX in NMOSD patients.
##plot MOG EDSS
forest.meta(nmosdedssmetagen,
sortvar = yi,
predict = TRUE,
print.tau2 = FALSE,
leftlabs = c("Study", "g", "SE"))
print(nmosdedssmetagen)## Review: NMOSDEDSSmeta
##
## Number of studies combined: k = 15
##
## MD 95%-CI z p-value
## Random effects model -1.1439 [-1.6257; -0.6620] -4.65 < 0.0001
##
## Quantifying heterogeneity:
## tau^2 = 0.2571 [0.0000; 1.5773]; tau = 0.5071 [0.0000; 1.2559]
## I^2 = 29.4% [0.0%; 62.0%]; H = 1.19 [1.00; 1.62]
##
## Test of heterogeneity:
## Q d.f. p-value
## 19.83 14 0.1357
##
## Details on meta-analytical method:
## - Inverse variance method
## - Restricted maximum-likelihood estimator for tau^2
## - Q-profile method for confidence interval of tau^2 and tauCombined analysis of EDSS
Here we see that there is no difference between EDSS scores of patients pre and post RTX between MOG and NMOSD (p=0.4254)
##plot mog and nmosd together
mergeedssmetagen <- update.meta(mergeedssmetagen,
subgroup = Group,
tau.common = FALSE)
forest.meta(mergeedssmetagen,
sortvar = yi,
predict = TRUE,
print.tau2 = FALSE,
subgroup = TRUE,
leftlabs = c("Study", "MD", "SE"))
print(mergeedssmetagen)## Review: mergeEDSSmeta
##
## Number of studies combined: k = 25
##
## MD 95%-CI z p-value
## Random effects model -0.9903 [-1.3590; -0.6217] -5.27 < 0.0001
##
## Quantifying heterogeneity:
## tau^2 = 0.3443 [0.0215; 1.2263]; tau = 0.5868 [0.1468; 1.1074]
## I^2 = 42.2% [6.9%; 64.1%]; H = 1.32 [1.04; 1.67]
##
## Test of heterogeneity:
## Q d.f. p-value
## 41.51 24 0.0146
##
## Results for subgroups (random effects model):
## k MD 95%-CI tau^2 tau Q I^2
## Group = MOG 10 -0.8384 [-1.4146; -0.2621] 0.4549 0.6744 20.20 55.4%
## Group = NMOSD 15 -1.1439 [-1.6257; -0.6620] 0.2571 0.5071 19.83 29.4%
##
## Test for subgroup differences (random effects model):
## Q d.f. p-value
## Between groups 0.64 1 0.4254
##
## Details on meta-analytical method:
## - Inverse variance method
## - Restricted maximum-likelihood estimator for tau^2
## - Q-profile method for confidence interval of tau^2 and tauMultiple regression - Univariate models
Here we removed mean baseline EDSS and prop children. As they significantly correlate with other predictors.
##meta-regression
mogreg<- read.csv("/Users/bosun/Documents/DPhil/OANG/6\ -\ Collaborations/2021-12-27_Val_meta/data/2022-03-27-mog_meta-regression data_last.csv")
mogreg <- mogreg %>%
dplyr::left_join(mogmeta, by = c("STUDY", "n"))
mogreg$weights <- 1/sqrt(mogreg$vi)
plot_metareg <- function(data = mogreg, variable = NULL){
mogreg2 <- mogreg[!is.na(mogreg[[paste0(variable)]]),]
reg <- rma(yi = yi,
vi = vi,
mods = mogreg2[[paste0(variable)]],
data = mogreg2)
# Specify basic plot, mapping sex to the x-axis, effect size 'd' to the y-axis,
# and 'weights' to the weight parameter.
preds <- predict.rma(reg)
plotdb <- cbind(mogreg2, preds)
ggplot(plotdb, aes(x = plotdb[[paste0(variable)]], y = plotdb[["yi"]], label=plotdb[["STUDY"]])) +
geom_ribbon(aes(ymin = ci.lb, ymax = ci.ub), alpha=0.05)+
geom_line(data = plotdb,aes(x = plotdb[[paste0(variable)]],y = plotdb$pred, size = 0.5), linetype="dashed", alpha=0.1)+
geom_text_repel(min.segment.length = 0, seed = 42, box.padding = 0.5, size=3, force=1)+
geom_point(aes(size=weights),shape = 21, alpha=1, fill="white") +
theme(legend.position = "none")+
xlab(variable)+
ylab("MD")
}
#Get modifiers
plotparam<- c("mean_time_to_RTX_months",
"mean_onset_age",
"mean_fup_months",
"mean_baseline_ARR",
"mean_baseline_EDSS",
"prop_relapsing",
"prop_children",
"prop_female")
#filter modifiers for highly correlated
mogreg[,unlist(plotparam)] %>%
chart.Correlation()
##remove prop_children as correlated with mean age onset
## remove mean EDSS as correlated highly with age of onset
plotparam <-plotparam[!plotparam %in% c("mean_baseline_EDSS", "prop_children")]
reg_plotlist <- lapply(c(plotparam,"n"), function(x){
p <- plot_metareg(data=mogreg, variable = x)
return(p)
})##plot bubble plots
patchwork::wrap_plots(reg_plotlist)
##get univariate model metrics
get_unimodels <- function(data, mod){
m_reg <- rma.uni(yi = yi,
vi = vi,
mods = data[[paste0(mod)]],
data = data)
metrics <- c(mod = mod, tau2 = round(m_reg$tau2,3),
tau =round(sqrt(m_reg$tau2), 3),
I2 = round(m_reg$I2, 3),
R2 = round(m_reg$R2, 3),
beta = round(m_reg$beta["mods",],3),
zval = round(m_reg$zval[2],3),
pval = round(m_reg$pval[2],7),
ci.lb = round(m_reg$ci.lb[2],3),
ci.ub = round(m_reg$ci.ub[2],3))
}
modregs <- lapply(c(plotparam,"n"), function(x){
get_unimodels(mogreg, x)
})
modregs <- arrange(as.data.frame(do.call(rbind, modregs)), desc(as.numeric(R2)))
kableExtra::kable(modregs)%>% kableExtra::kable_material(c("striped", "hover"))| mod | tau2 | tau | I2 | R2 | beta.mods | zval | pval | ci.lb | ci.ub |
|---|---|---|---|---|---|---|---|---|---|
| mean_time_to_RTX_months | 0 | 0 | 0 | 100 | 0.016 | 4.118 | 3.82e-05 | 0.008 | 0.024 |
| mean_baseline_ARR | 0 | 0 | 0 | 100 | -0.85 | -5.035 | 5e-07 | -1.18 | -0.519 |
| prop_female | 0.196 | 0.443 | 53.192 | 50.714 | -3.67 | -1.922 | 0.0546228 | -7.413 | 0.073 |
| n | 0.194 | 0.441 | 48.954 | 10.659 | -0.006 | -0.954 | 0.3402878 | -0.017 | 0.006 |
| mean_onset_age | 0.303 | 0.55 | 59.152 | 0 | 0.028 | 1.519 | 0.1286385 | -0.008 | 0.063 |
| mean_fup_months | 0.239 | 0.489 | 59.02 | 0 | -0.008 | -0.651 | 0.5151981 | -0.03 | 0.015 |
| prop_relapsing | 0.244 | 0.494 | 58.782 | 0 | -0.978 | -0.692 | 0.4889695 | -3.747 | 1.791 |
Multi-model inference
Here, all possible predictor combinations with multi-model inference, we see that mean ARR is the most important predictor of RTX response.
multimodel.inference(TE = "yi",
seTE = "sei",
data = mogreg,
predictors = c(plotparam),
interaction = TRUE)## You entered 6 predictors. Interactions can only be modeled for four or less predictors. Therefore, no interactions are modeled.
##
|
| | 0%
|
|= | 2%
|
|== | 3%
|
|=== | 5%
|
|==== | 6%
|
|===== | 8%
|
|======= | 9%
|
|======== | 11%
|
|========= | 12%
|
|========== | 14%
|
|=========== | 16%
|
|============ | 17%
|
|============= | 19%
|
|============== | 20%
|
|=============== | 22%
|
|================ | 23%
|
|================== | 25%
|
|=================== | 27%
|
|==================== | 28%
|
|===================== | 30%
|
|====================== | 31%
|
|======================= | 33%
|
|======================== | 34%
|
|========================= | 36%
|
|========================== | 38%
|
|=========================== | 39%
|
|============================ | 41%
|
|============================== | 42%
|
|=============================== | 44%
|
|================================ | 45%
|
|================================= | 47%
|
|================================== | 48%
|
|=================================== | 50%
|
|==================================== | 52%
|
|===================================== | 53%
|
|====================================== | 55%
|
|======================================= | 56%
|
|======================================== | 58%
|
|========================================== | 59%
|
|=========================================== | 61%
|
|============================================ | 62%
|
|============================================= | 64%
|
|============================================== | 66%
|
|=============================================== | 67%
|
|================================================ | 69%
|
|================================================= | 70%
|
|================================================== | 72%
|
|=================================================== | 73%
|
|==================================================== | 75%
|
|====================================================== | 77%
|
|======================================================= | 78%
|
|======================================================== | 80%
|
|========================================================= | 81%
|
|========================================================== | 83%
|
|=========================================================== | 84%
|
|============================================================ | 86%
|
|============================================================= | 88%
|
|============================================================== | 89%
|
|=============================================================== | 91%
|
|================================================================= | 92%
|
|================================================================== | 94%
|
|=================================================================== | 95%
|
|==================================================================== | 97%
|
|===================================================================== | 98%##
##
## Multimodel Inference: Final Results
## --------------------------
##
## - Number of fitted models: 64
## - Full formula: ~ mean_time_to_RTX_months + mean_onset_age + mean_fup_months + mean_baseline_ARR + prop_relapsing + prop_female
## - Coefficient significance test: knha
## - Interactions modeled: no
## - Evaluation criterion: AICc
##
##
## Best 5 Models
## --------------------------
##
##
## Global model call: metafor::rma(yi = TE, sei = seTE, mods = form, data = glm.data,
## method = method, test = test)
## ---
## Model selection table
## (Int) men_bsl_ARR men_fup_mnt men_ons_age prp_fml prp_rlp df logLik AICc
## 2 + -0.8497 3 -5.647 19.7
## 18 + -0.7958 -1.907 4 -2.903 21.8
## 4 + -0.8756 0.006311 4 -4.954 22.4
## 6 + -0.8607 0.01268 4 -4.452 22.6
## 34 + -0.8387 -0.7559 4 -5.177 22.8
## delta weight
## 2 0.00 0.486
## 18 2.11 0.169
## 4 2.66 0.129
## 6 2.92 0.113
## 34 3.11 0.103
## Models ranked by AICc(x)
##
##
## Multimodel Inference Coefficients
## --------------------------
##
##
## Estimate Std. Error z value Pr(>|z|)
## intrcpt 0.3983073530 0.602707328 0.66086363 0.5086998
## mean_baseline_ARR -0.8381088054 0.150244981 5.57828156 0.0000000
## prop_female -0.3178070621 0.800662121 0.39693031 0.6914189
## mean_fup_months 0.0008785040 0.002877126 0.30534082 0.7601066
## mean_onset_age 0.0016179293 0.005125106 0.31568701 0.7522401
## prop_relapsing -0.0923513767 0.348598158 0.26492216 0.7910694
## mean_time_to_RTX_months 0.0001573002 0.001721205 0.09138954 0.9271831
##
##
## Predictor Importance
## --------------------------
##
##
## model importance
## 1 mean_baseline_ARR 0.98460599
## 2 prop_female 0.16270696
## 3 mean_fup_months 0.14308581
## 4 mean_onset_age 0.12012746
## 5 prop_relapsing 0.11775970
## 6 mean_time_to_RTX_months 0.08825281
Subgroup extraction and analysis
## $relapsing_mog
## Review: relapsing_mog
##
## Number of studies combined: k = 16
##
## MD 95%-CI z p-value
## Random effects model -1.0270 [-1.3226; -0.7314] -6.81 < 0.0001
##
## Quantifying heterogeneity:
## tau^2 = 0.1369 [0.0045; 0.8841]; tau = 0.3700 [0.0673; 0.9403]
## I^2 = 47.6% [6.3%; 70.7%]; H = 1.38 [1.03; 1.85]
##
## Test of heterogeneity:
## Q d.f. p-value
## 28.60 15 0.0181
##
## Details on meta-analytical method:
## - Inverse variance method
## - Restricted maximum-likelihood estimator for tau^2
## - Q-profile method for confidence interval of tau^2 and tau
##
## $adult_mog
## Review: adult_mog
##
## Number of studies combined: k = 10
##
## MD 95%-CI z p-value
## Random effects model -1.0018 [-1.5620; -0.4416] -3.50 0.0005
##
## Quantifying heterogeneity:
## tau^2 = 0.4621 [0.0956; 1.4413]; tau = 0.6798 [0.3092; 1.2006]
## I^2 = 85.5% [75.2%; 91.6%]; H = 2.63 [2.01; 3.44]
##
## Test of heterogeneity:
## Q d.f. p-value
## 62.27 9 < 0.0001
##
## Details on meta-analytical method:
## - Inverse variance method
## - Restricted maximum-likelihood estimator for tau^2
## - Q-profile method for confidence interval of tau^2 and tau
##
## $children_mog
## Review: children_mog
##
## Number of studies combined: k = 7
##
## MD 95%-CI z p-value
## Random effects model -1.1957 [-1.8592; -0.5322] -3.53 0.0004
##
## Quantifying heterogeneity:
## tau^2 = 0.5047 [0.0592; 6.4154]; tau = 0.7105 [0.2432; 2.5329]
## I^2 = 65.2% [21.7%; 84.5%]; H = 1.69 [1.13; 2.54]
##
## Test of heterogeneity:
## Q d.f. p-value
## 17.23 6 0.0085
##
## Details on meta-analytical method:
## - Inverse variance method
## - Restricted maximum-likelihood estimator for tau^2
## - Q-profile method for confidence interval of tau^2 and tau
##
## $fu12_mog
## Review: fu12_mog
##
## Number of studies combined: k = 14
##
## MD 95%-CI z p-value
## Random effects model -1.1253 [-1.5243; -0.7264] -5.53 < 0.0001
##
## Quantifying heterogeneity:
## tau^2 = 0.2533 [0.0100; 1.2335]; tau = 0.5033 [0.0998; 1.1106]
## I^2 = 50.3% [8.2%; 73.1%]; H = 1.42 [1.04; 1.93]
##
## Test of heterogeneity:
## Q d.f. p-value
## 26.15 13 0.0162
##
## Details on meta-analytical method:
## - Inverse variance method
## - Restricted maximum-likelihood estimator for tau^2
## - Q-profile method for confidence interval of tau^2 and tau
##
## $first_mog
## Review: firstline_mog
##
## Number of studies combined: k = 5
##
## MD 95%-CI z p-value
## Random effects model -1.3202 [-2.6540; 0.0136] -1.94 0.0524
##
## Quantifying heterogeneity:
## tau^2 = 1.8373 [0.4004; 17.6857]; tau = 1.3555 [0.6327; 4.2054]
## I^2 = 84.9% [66.3%; 93.2%]; H = 2.57 [1.72; 3.84]
##
## Test of heterogeneity:
## Q d.f. p-value
## 26.44 4 < 0.0001
##
## Details on meta-analytical method:
## - Inverse variance method
## - Restricted maximum-likelihood estimator for tau^2
## - Q-profile method for confidence interval of tau^2 and tauForest plots of subgroups
for (i in 1:length(subgroupsmetagen)){
forest.meta(subgroupsmetagen[[paste0(names(subgroupsmetagen[i]))]],
sortvar = yi,
predict = TRUE,
print.tau2 = FALSE,
leftlabs = c("Study", "g", "SE"))
grid::grid.text(paste0(names(subgroupsmetagen[i])), .5, .9, gp=grid::gpar(cex=2))
}
